Lecture 1
Introduction
Outline
- What are INLA and
inlabru? - Why the Bayesian framework?
- Which model are
inlabru-friendly? - What are Latent Gaussian Models?
- How are they implemented in
inlabru?
What is INLA? What is inlabru?
The short answer:
INLA is a fast method to do Bayesian inference with latent Gaussian models and
inlabruis anR-package that implements this method with a flexible and simple interface.
The (much) longer answer:
- Rue, H., Martino, S. and Chopin, N. (2009), Approximate Bayesian inference for latent Gaussian models by using integrated nested Laplace approximations. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 71: 319-392.
- Van Niekerk, J., Krainski, E., Rustand, D., & Rue, H. (2023). A new avenue for Bayesian inference with INLA. Computational Statistics & Data Analysis, 181, 107692.
- Lindgren, F., Bachl, F., Illian, J., Suen, M. H., Rue, H., & Seaton, A. E. (2024). inlabru: software for fitting latent Gaussian models with non-linear predictors. arXiv preprint arXiv:2407.00791.
- Lindgren, F., Bolin, D., & Rue, H. (2022). The SPDE approach for Gaussian and non-Gaussian fields: 10 years and still running. Spatial Statistics, 50, 100599.
Where?
So… Why should you use inlabru?
- What type of problems can we solve?
- What type of models can we use?
- When can we use it?
So… Why should you use inlabru?
- What type of problems can we solve?
- What type of models can we use?
- When can we use it `inlabru’?
To give proper answers to these questions, we need to start at the very beginning…
The core
- We have observed something.
The core
- We have observed something.
- We have questions.
The core
- We have observed something.
- We have questions.
- We want answers!
How do we find answers?
We need to make choices:
- Bayesian or frequentist?
- How do we model the data?
- How do we compute the answer?
How do we find answers?
We need to make choices:
- Bayesian or frequentist?
- How do we model the data?
- How do we compute the answer?
These questions are not independent.
Bayesian or frequentist?
In this course we embrace the Bayesian perspective
- There are no “true but unknown” parameters !
Bayesian or frequentist?
In this course we embrace the Bayesian perspective
- There are no “true but unknown” parameters!
- Every parameter is described by a probability distribution!
Bayesian or frequentist?
In this course we embrace the Bayesian perspective
- There are no “true but unknown” parameters!
- Every parameter is described by a probability distribution!
- Evidence from the data is used to update the belief we had before observing the data!
Some more details I
We define as the linear predictor the mean (or a function of the mean) of our observations given the model components.
In this case \(E(y_i|\beta_0, \beta_i) = \eta_i = \beta_0 + \beta_1 x_i\)
Some more details I
We define as the linear predictor the mean (or a function of the mean) of our observations given the model components.
In this case \(E(y_i|\beta_0, \beta_i) =\eta_i = \color{red}{\boxed{\beta_0}} + \color{red}{\boxed{\beta_1 x_i}}\)
This model has two components!
Some more details II
Given the linear predictor \(\eta\) the observations are independent of each other!
This means that all dependencies in the observations are accounted for by the components!
Some more details II
Given the linear predictor \(\eta\) the observations are independent of each other!
The observation model (likelihood) can be written as:
\[ \pi(\mathbf{y}|\eta,\sigma^2) = \prod_{i = 1}^n\pi(y_i|\eta_i,\sigma^2) \]
Let’s formalize this a bit…
The elements of an inlabru-friendly statistical model are:
- The observational model \[ \begin{aligned} y_i|\eta_i, \sigma^2 & \sim\mathcal{N}(\eta_i,\sigma^2),\qquad i = 1,\dots,n\\ E(y_i|\eta_i, \sigma^2) & = \eta_i \end{aligned} \] Note: We assume that, given the linear predictor \(\eta\), the data are independent of each other! Data dependence is expressed through the components of the linear predictor.
Let’s formalize this a bit…
The elements of an inlabru-friendly statistical model are:
The observational model \(y_i|\eta_i,\sigma^2\sim\mathcal{N}(\eta_i,\sigma^2),\qquad i = 1,\dots,n\)
A model for the linear predictor \[ E(y_i|\eta_i,\sigma^2) = \eta_i = \beta_0 + \beta_1x_i \]
Let’s formalize this a bit…
The elements of a inlabru friendly statistical model are:
The observational model \(y_i|\eta_i,\sigma^2\sim\mathcal{N}(\eta_i,\sigma^2),\qquad i = 1,\dots,n\)
A model for the linear predictor
\[ E(y_i|\eta_i,\sigma^2) = \eta_i = \color{red}{\boxed{\beta_0}} + \color{red}{\boxed{\beta_1x_i} } \]
Note 1: These are the components of our model! These explain the dependence structure of the data.
Let’s formalize this a bit…
The elements of a inlabru friendly statistical model are:
The observational model \(y_i|\eta_i,\sigma^2\sim\mathcal{N}(\eta_i,\sigma^2),\qquad i = 1,\dots,n\)
A model for the linear predictor \(\eta_i = \color{red}{\boxed{\beta_0}} + \color{red}{\boxed{\beta_1x_i} }\)
A prior for the model components \(\textbf{u}\) \[ \mathbf{u} = \{\beta_0, \beta_1\}\sim\mathcal{N}(0,\mathbf{Q}^{-1}) \] Note: These always have a Gaussian prior and are used to explain the dependencies among observations!
Let’s formalize this a bit…
The elements of a inlabru friendly statistical model are:
The observational model \(y_i|\eta_i,\sigma^2\sim\mathcal{N}(\eta_i,\sigma^2),\qquad i = 1,\dots,n\)
A model for the linear predictor \(\eta_i = \color{red}{\boxed{\beta_0}} + \color{red}{\boxed{\beta_1x_i} }\)
A prior for the model components \(\mathbf{u} = \{\beta_0, \beta_1\}\sim\mathcal{N}(0,\mathbf{Q}^{-1})\)
A prior for the non-Gaussian parameters \(\theta\) \[ \theta = \sigma^2 \]
Latent Gaussian Models (LGM)
The observation model: \[ \pi(\mathbf{y}|\eta,\theta) = \prod_{i=1}^{n}\pi(y_i|\eta_i,\theta) \]
Linear predictor \(\eta_i = \beta_0 + \beta_1 x_i\)
Latent Gaussian field \(\pi(\mathbf{u}|\theta)\)
The hyperparameters: \(\pi(\theta)\)
- Stage 1 The data generating process
Latent Gaussian Models (LGM)
The observation model: \[ \pi(\mathbf{y}|\eta,\theta) = \prod_{i=1}^{n}\pi(y_i|\eta_i,\theta) \]
Linear predictor \(\eta_i = \beta_0 + \beta_1 x_i\)
Latent Gaussian field \(\pi(\mathbf{u}|\theta)\)
The hyperparameters: \(\pi(\theta)\)
Stage 1 The data generating process
Stage 2 The dependence structure
Latent Gaussian Models (LGM)
The observation model: \[ \pi(\mathbf{y}|\eta,\theta) = \prod_{i=1}^{n}\pi(y_i|\eta_i,\theta) \]
Linear predictor \(\eta_i = \beta_0 + \beta_1 x_i\)
Latent Gaussian field \(\pi(\mathbf{u}|\theta)\)
The hyperparameters: \(\pi(\theta)\)
Stage 1 The data generating process
Stage 2 The dependence structure
Stage 3 The hyperparameters
Latent Gaussian Models (LGM)
The observation model: \[ \pi(\mathbf{y}|\eta,\theta) = \prod_{i=1}^{n}\pi(y_i|\eta_i,\theta) \]
Linear predictor \(\eta_i = \beta_0 + \beta_1 x_i\)
Latent Gaussian field \(\pi(\mathbf{u}|\theta)\)
The hyperparameters: \(\pi(\theta)\)
Stage 1 The data generating process
Stage 2 The dependence structure
Stage 3 The hyperparameters
Q: What are we interested in?
The posterior distribution
\[ \color{red}{\pi(\mathbf{u},\theta|\mathbf{y})}\propto \color{blue}{\pi(\mathbf{y}|\mathbf{u},\theta)}\color{green}{\pi(\mathbf{u}|\theta)\pi(\theta)} \]
The posterior distribution
inlabru for linear regression
The Model \[ \begin{aligned} y_i|\eta_i, \sigma^2 & \sim \mathcal{N}(\eta_i,\sigma^2)\\ \eta_i & = \beta_0 + \beta_i x_i \end{aligned} \]
library(inlabru)
Pallid[1:3,c("w","tl")] w tl
1 2.239 95.9
2 2.948 95.0
3 3.402 108.0
The code
# define model components
cmp = ~ -1 + beta0(1) + beta1(tl, model = "linear")
# define model predictor
eta = w ~ beta0 + beta1
# build the observation model
lik = bru_obs(formula = eta,
family = "gaussian",
data = Pallid)
# fit the model
fit = bru(cmp, lik)inlabru for linear regression
The Model \[ \begin{aligned} y_i|\eta_i, \sigma^2 & \sim \mathcal{N}(\eta_i,\sigma^2)\\ \eta_i & = \color{red}{\boxed{\beta_0}} + \color{red}{\boxed{\beta_i x_i}} \end{aligned} \]
Pallid[1:3,c("w","tl")] w tl
1 2.239 95.9
2 2.948 95.0
3 3.402 108.0
The code
inlabru for linear regression
The Model \[ \begin{aligned} y_i|\eta_i, \sigma^2 & \sim \mathcal{N}(\eta_i,\sigma^2)\\ \eta_i & = \color{red}{\boxed{\beta_0 + \beta_i x_i}} \end{aligned} \]
library(inlabru)
Pallid[1:3,c("w","tl")] w tl
1 2.239 95.9
2 2.948 95.0
3 3.402 108.0
The code
inlabru for linear regression
The Model \[ \begin{aligned} \color{red}{\boxed{y_i|\eta_i, \sigma^2}} & \color{red}{\boxed{\sim \mathcal{N}(\eta_i,\sigma^2)}}\\ \eta_i & = \beta_0 + \beta_i x_i \end{aligned} \]
library(inlabru)
Pallid[1:3,c("w","tl")] w tl
1 2.239 95.9
2 2.948 95.0
3 3.402 108.0
The code
inlabru for linear regression
The Model \[ \begin{aligned} y_i|\eta_i, \sigma^2 & \sim \mathcal{N}(\eta_i,\sigma^2)\\ \eta_i & = \beta_0 + \beta_i x_i \end{aligned} \]
library(inlabru)
Pallid[1:3,c("w","tl")] w tl
1 2.239 95.9
2 2.948 95.0
3 3.402 108.0
The code
inlabru for linear regression
Real datasets are more complicated!
Data can have several dependence structures: temporal, spatial,…
Using a Bayesian framework:
Build (hierarchical) models to account for potentially complicated dependency structures in the data.
Attribute uncertainty to model parameters and latent variables using priors.
Two main challenges:
- Need computationally efficient methods to calculate posteriors (this is where INLA helps!).
- Select priors in a sensible way (we’ll talk about this)
The good news!!
In many cases complicated spatio-temporal models are just special cases of the same model structure!! 😃
Stage 1: What is the distribution of the responses?
Stage 2: What are the model components? and what is their distribution?
Stage 3: What are our prior beliefs about the parameters controlling the components in the model?
The good news!!
In many cases complicated spatio-temporal models are just special cases of the same model structure!! 😃
Stage 1: What is the distribution of the responses?
- Gaussian response? (temperature, rainfall, fish weight …)
- Count data? (people infected with a disease in each area)
- Point pattern? (locations of trees in a forest)
- Binary data? (yes/no response, binary image)
- Survival data? (recovery time, time to death)
- … (many more examples!!)
Stage 2: What are the model components? and what is their distribution?
Stage 3: What are our prior beliefs about the parameters controlling the components in the model?
The good news!!
In many cases complicated spatio-temporal models are just special cases of the same model structure!! 😃
Stage 1: What is the distribution of the responses?
- We assume data to be conditionally independent given the model components and some hyperparameters
- This means that all dependencies in data are explained in Stage
Stage 2: What are the model components? and what is their distribution?
Stage 3: What are our prior beliefs about the parameters controlling the components in the model?
The good news!!
In many cases complicated spatio-temporal models are just special cases of the same model structure!! 😃
Stage 1: What is the distribution of the responses?
Stage 2: What are the model components? and what is their distribution?
Here we can have:
- Fixed effects for covariates
- Unstructured random effects (individual effects, group effects)
- Structured random effects (AR(1), regional effects, )
- …
These are linked to the responses in the likelihood through linear predictors.
- Stage 3: What are our prior beliefs about the parameters controlling the components in the model?
The good news!!
In many cases complicated spatio-temporal models are just special cases of the same model structure!! 😃
Stage 1: What is the distribution of the responses?
Stage 2: What are the model components? and what is their distribution?
Stage 3: What are our prior beliefs about the parameters controlling the components in the model?
The likelihood and the latent model typically have hyperparameters that control their behavior.
They can include:
- Variance of observation noise
- Dispersion parameter in the negative binomial model
- Variance of unstructured effects
- …
The second good news!
No matter how complicated your model is, the inlabru workflow is always the same 😃
The second good news!
No matter how complicated your model is, the inlabru workflow is always the same 😃
NOTE we will see later that this function can also be non-linear….😁
Examples
The Tokyo rainfall data
One example with time series: Rainfall over 1 mm in the Tokyo area for each calendar day during two years (1983-84)
The model
Stage 1 The observation model
\[ y_t|\eta_t\sim\text{Bin}(n_t, p_t),\qquad \eta_t = \text{logit}(p_t),\qquad i = 1,\dots,366 \] \[ n_t = \left\{ \begin{array}{lr} 1, & \text{for}\; 29\; \text{February}\\ 2, & \text{other days} \end{array}\right. \] \[ y_t = \begin{cases} \{0,1\}, & \text{for}\; 29\; \text{February}\\ \{0,1,2\}, & \text{other days} \end{cases} \]
- the likelihood has no hyperparameters
The model
Stage 1 The observation model
\[ y_t|\eta_t\sim\text{Bin}(n_t, p_t),\qquad \eta_t = \text{logit}(p_t),\qquad i = 1,\dots,366 \]
Stage 2 The latent field \[ \eta_t = \beta_0 + f(\text{time}_t) \]
probability of rain depends on the day of the year \(t\)
\(\beta_0\) is an intercept
\(f(\text{time}_t)\) is a RW2 model (this is just a smoother). The smoothness is controlled by a hyperparameter \(\tau_f\)
The model
Stage 1 The observation model
\[ y_t|\eta_t\sim\text{Bin}(n_t, p_t),\qquad \eta_t = \text{logit}(p_t),\qquad i = 1,\dots,366 \]
Stage 2 The latent field \[ \eta_t = \beta_0 + f(\text{time}_t) \]
Stage 3 The hyperparameters
The structured time effect is controlled by one parameter \(\tau_f\).
We assign a prior to \(\tau_f\) to finalize the model.
inlabru for time series
The Model
\[ \begin{aligned} y_t|\eta_t & \sim \text{Binomial}(n_t,p_t)\\ \text{logit}(p_t) = \eta_i & = \color{red}{\boxed{\beta_0}} + \color{red}{\boxed{f(\text{time}_t)}} \end{aligned} \]
Tokyo[1:3,] y n time
1 0 2 1
2 0 2 2
3 1 2 3
The code
inlabru for time series
The Model
\[ \begin{aligned} y_t|\eta_t & \sim \text{Binomial}(n_t,p_t)\\ \text{logit}(p_t) = \color{red}{\boxed{\eta_i}} & = \color{red}{\boxed{\beta_0 + f(\text{time}_t)}} \end{aligned} \]
Tokyo[1:3,] y n time
1 0 2 1
2 0 2 2
3 1 2 3
The code
inlabru for time series
The Model
\[ \begin{aligned} \color{red}{\boxed{y_t|\eta_t}} & \color{red}{\boxed{\sim \text{Binomial}(n_t,p_t)}}\\ \text{logit}(p_t) = \eta_i & = \beta_0 + f(\text{time}_t) \end{aligned} \]
Tokyo[1:3,] y n time
1 0 2 1
2 0 2 2
3 1 2 3
The code
inlabru for time series
The Model
\[ \begin{aligned} y_t|\eta_t & \sim \text{Binomial}(n_t,p_t)\\ \text{logit}(p_t) = \eta_i & = \beta_0 + f(\text{time}_t) \end{aligned} \]
Tokyo[1:3,] y n time
1 0 2 1
2 0 2 2
3 1 2 3
The code
Example: disease mapping
We observed larynx cancer mortality counts for males in 544 district of Germany from 1986 to 1990 and want to make a model.
\(y_i\): The count at location \(i\).
\(E_i\): An offset; expected number of cases in district \(i\).
\(c_i\): A covariate (level of smoking consumption) at \(i\)
\(\boldsymbol{s}_i\): spatial location \(i\) .
Bayesian disease mapping
Stage 1: We assume the responses are Poisson distributed: \[ y_i \mid \eta_i \sim \text{Poisson}(E_i\exp(\eta_i))) \]
Stage 2: \(\eta_i\) is a linear function of three components: an intercept, a covariate \(c_i\), a spatially structured effect \(\omega\) likelihood by \[ \eta_i = \beta_0 + \beta_1\ c_i + \omega_i \]
Stage 3:
- \(\tau_{\omega}\): Precisions parameter for the random effects
. . .
The latent field is \(\boldsymbol{u} = (\beta_0, \beta_1, \omega_1, \omega_2,\ldots, \omega_n)\), the hyperparameters are \(\boldsymbol{\theta} = (\tau_{\omega})\), and must be given a prior.
inlabru for disease mapping
The Model
\[ \begin{aligned} y_i|\eta_t & \sim \text{Poisson}(E_i\lambda_i)\\ \text{log}(\lambda_i) = \eta_i & = \color{red}{\boxed{\beta_0}} + \color{red}{\boxed{\beta_1\ c_i}} + \color{red}{\boxed{\omega_i}} \end{aligned} \]
g = system.file("demodata/germany.graph",
package="INLA")
Germany[1:3,] region E Y x region.struct
1 1 7.965008 8 56 1
2 2 22.836219 22 65 2
3 3 22.094716 19 50 3
The code
# define model component
cmp = ~ -1 + beta0(1) + beta1(x, model = "linear") +
space(region, model = "besag", graph = g)
# define model predictor
eta = Y ~ beta0 + beta1 + space
# build the observation model
lik = bru_obs(formula = eta,
family = "poisson",
E = E,
data = Germany)
# fit the model
fit = bru(cmp, lik)inlabru for disease mapping
The Model
\[ \begin{aligned} y_i|\eta_t & \sim \text{Poisson}(E_i\lambda_i)\\ \text{log}(\lambda_i) = \color{red}{\boxed{\eta_i}} & = \color{red}{\boxed{\beta_0 + \beta_1\ c_i + \omega_i}} \end{aligned} \]
g = system.file("demodata/germany.graph",
package="INLA")
Germany[1:3,] region E Y x region.struct
1 1 7.965008 8 56 1
2 2 22.836219 22 65 2
3 3 22.094716 19 50 3
The code
# define model component
cmp = ~ -1 + beta0(1) + beta1(x, model = "linear") +
space(region, model = "bym2", graph = g)
# define model predictor
eta = Y ~ beta0 + beta1 + space
# build the observation model
lik = bru_obs(formula = eta,
family = "poisson",
E = E,
data = Germany)
# fit the model
fit = bru(cmp, lik)inlabru for disease mapping
The Model
\[ \begin{aligned} \color{red}{\boxed{y_i|\eta_t}} & \color{red}{\boxed{\sim \text{Poisson}(E_i\lambda_i)}}\\ \text{log}(\lambda_i) = \eta_i & = \beta_0 + \beta_1\ c_i + \omega_i \end{aligned} \]
g = system.file("demodata/germany.graph",
package="INLA")
Germany[1:3,] region E Y x region.struct
1 1 7.965008 8 56 1
2 2 22.836219 22 65 2
3 3 22.094716 19 50 3
The code
# define model component
cmp = ~ -1 + beta0(1) + beta1(x, model = "linear") +
space(region, model = "bym2", graph = g)
# define model predictor
eta = Y ~ beta0 + beta1 + space
# build the observation model
lik = bru_obs(formula = eta,
family = "poisson",
E = E,
data = Germany)
# fit the model
fit = bru(cmp, lik)inlabru for disease mapping
Bayesian Geostatistics
Encounter probability of Pacific Cod (Gadus macrocephalus) from a trawl survey.
\(y(s)\) Presence or absence in location \(s\)
Bayesian Geostatistics
Stage 1 Model for the response \[ y(s)|\eta(s)\sim\text{Binom}(1, p(s)) \]
Stage 2 Latent field model \[ \eta(s) = \text{logit}(p(s)) = \beta_0 + f( x(s)) + \omega(s) \]
Stage 3 Hyperparameters
Bayesian Geostatistics
- Stage 1 Model for the response \[ y(s)|\eta(s)\sim\text{Binom}(1, p(s)) \]
- Stage 2 Latent field model \[
\eta(s) = \text{logit}(p(s)) = \beta_0 + f( x(s)) + \omega(s)
\]
- A global intercept \(\beta_0\)
- A smooth effect of covariate \(x(s)\) (depth)
- A Gaussian field \(\omega(s)\) (will discuss this later..)
- Stage 3 Hyperparameters
Bayesian Geostatistics
Stage 1 Model for the response \[ y(s)|\eta(s)\sim\text{Binom}(1, p(s)) \]
Stage 2 Latent field model \[ \eta(s) = \text{logit}(p(s)) = \beta_0 + \beta_1 x(s) + \omega(s) \]
Stage 3 Hyperparameters
- Precision for the smooth function \(f(\cdot)\)
- Range and sd in the Gaussian field \(\sigma_{\omega}, \tau_{\omega}\)
inlabru for geostatistics
The Model
\[ \begin{aligned} y(s)|\eta(s) & \sim\text{Binom}(1, p(s))\\ \eta(s) & = \color{red}{\boxed{\beta_0}} + \color{red}{\boxed{ f(x(s))}} + \color{red}{\boxed{ \omega(s)}}\\ \end{aligned} \]
df %>% select(depth, present) %>% print(n = 3)Simple feature collection with 2143 features and 2 fields
Geometry type: POINT
Dimension: XY
Bounding box: xmin: 343.0617 ymin: 5635.893 xmax: 579.3681 ymax: 5839.019
Projected CRS: +proj=utm +zone=9 +datum=WGS84 +no_defs +type=crs +units=km
# A tibble: 2,143 × 3
depth present geometry
<dbl> <dbl> <POINT [km]>
1 201 1 (446.4752 5793.426)
2 212 1 (446.4594 5800.136)
3 220 0 (448.5987 5801.687)
# ℹ 2,140 more rows
The code
# define model component
cmp = ~ -1 + Intercept(1) + depth_smooth(log(depth), model='rw2') +
space(geometry, model = spde_model)
# define model predictor
eta = present ~ Intercept + depth_smooth + space
# build the observation model
lik = bru_obs(formula = eta,
data = df,
family = "binomial")
# fit the model
fit = bru(cmp, lik)inlabru for geostatistics
The Model
\[ \begin{aligned} y(s)|\eta(s) & \sim\text{Binom}(1, p(s))\\ \color{red}{\boxed{\eta(s)}} & = \color{red}{\boxed{\beta_0 + f(x(s)) + \omega(s)}}\\ \end{aligned} \]
df %>% select(depth, present) %>% print(n = 3)Simple feature collection with 2143 features and 2 fields
Geometry type: POINT
Dimension: XY
Bounding box: xmin: 343.0617 ymin: 5635.893 xmax: 579.3681 ymax: 5839.019
Projected CRS: +proj=utm +zone=9 +datum=WGS84 +no_defs +type=crs +units=km
# A tibble: 2,143 × 3
depth present geometry
<dbl> <dbl> <POINT [km]>
1 201 1 (446.4752 5793.426)
2 212 1 (446.4594 5800.136)
3 220 0 (448.5987 5801.687)
# ℹ 2,140 more rows
The code
# define model component
cmp = ~ -1 + Intercept(1) + depth_smooth(log(depth), model='rw2') +
space(geometry, model = spde_model)
# define model predictor
eta = present ~ Intercept + depth_smooth + space
# build the observation model
lik = bru_obs(formula = eta,
data = df,
family = "binomial")
# fit the model
fit = bru(cmp, lik)inlabru for geostatistics
The Model
\[ \begin{aligned} \color{red}{\boxed{y(s)|\eta(s)}} & \sim \color{red}{\boxed{\text{Binom}(1, p(s))}}\\ \eta(s) & = \beta_0 + f(x(s)) + \omega(s)\\ \end{aligned} \]
df %>% select(depth, present) %>% print(n = 3)Simple feature collection with 2143 features and 2 fields
Geometry type: POINT
Dimension: XY
Bounding box: xmin: 343.0617 ymin: 5635.893 xmax: 579.3681 ymax: 5839.019
Projected CRS: +proj=utm +zone=9 +datum=WGS84 +no_defs +type=crs +units=km
# A tibble: 2,143 × 3
depth present geometry
<dbl> <dbl> <POINT [km]>
1 201 1 (446.4752 5793.426)
2 212 1 (446.4594 5800.136)
3 220 0 (448.5987 5801.687)
# ℹ 2,140 more rows
The code
# define model component
cmp = ~ -1 + Intercept(1) + depth_smooth(log(depth), model='rw2') +
space(geometry, model = spde_model)
# define model predictor
eta = present ~ Intercept + depth_smooth + space
# build the observation model
lik = bru_obs(formula = eta,
data = df,
family = "binomial")
# fit the model
fit = bru(cmp, lik)inlabru for geostatistics
Take home message!
Many of the models you have used (and some you have never used but will learn about) are just special cases of the large class of Latent Gaussian models
inlabru provides an efficient and unified way to fit all these models!